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ABSTRACT 

Spectral retrieval has proven to be a powerful tool for constraining the physical properties and atmo¬ 
spheric compositions of extrasolar planet atmospheres from observed spectra, primarily for transiting 
objects but also for directly imaged planets and brown dwarfs. Despite its strengths, this approach 
has been applied to only about a dozen targets. Determining the abundances of the main carbon 
and oxygen-bearing compounds in a planetary atmosphere can lead to the C/0 ratio of the object, 
which is crucial in understanding its formation and migration history. We present a retrieval analysis 
on the published near-infrared spectrum of k And b, a directly imaged substellar companion to a 
young B9 star. We fit the emission spectrum model utilizing a Markov Chain Monte Carlo algorithm. 

We estimate the abundance of water vapor, and its uncertainty, in the atmosphere of the object. In 
addition, we place an upper limit on the abundance of CH 4 . We compare qualitatively our results to 
studies that have applied model retrieval on multiband photometry and emission spectroscopy of hot 
Jupiters (extrasolar giant planets with orbital periods of several days) and the directly imaged giant 
planet HR 8799b. 

Subject headings: stars: planetary systems — direct imaging - techniques: spectroscopic — methods: 
data analysis — planets and satellites: individual (k And b) 


1. INTRODUCTION 

The study of exoplanets is moving from an era focussed 
on discovery to one of characterization. While most ex¬ 
oplanets detected to date have been discovered through 
the radial velocity or transit techniques, there are a hand¬ 
ful of substellar companions within ~ 100 AU of their 
host stars discovered through direct imaging (Pepe et al. 
2014). It remains unclear whether these populations are 
distinct or can be considered as part of a continuum of 
objects drawn from a variety of formation mechanisms, 
and experiencing a range of evolutionary histories. 

The compositions of transiting and directly imaged ex¬ 
trasolar giant planets may be an important marker of 
their formation and migration history (Oberg et al. 2011; 
Madhusudhan et al. 2014a). These authors suggest that 
giant planets’ C/0 and C/H ratios are affected by the 
locations in the protoplanetary disks where they formed 
with respect to the water and CO ice lines, as well as 
by the formation mechanism. E.g., Oberg et al. (2011) 
suggest that stellar C/0 and C/H ratios are indicative of 
either formation within the water ice line or via gravita¬ 
tional instability. Substellar C/0, but superstellar C/H 
point to significant accretion of water ice bodies after 
the gas accretion phase. Superstellar C/0 and C/H can 
mean gas accumulation near the CO or CO 2 ice lines or 
significant accretion of carbon grains. Superstellar C/0 
and substellar C/H indicates gas accretion from outside 
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the water snow line. A caveat is that Oberg et al. (2011) 
assume that the planet forms at a fixed location within 
the disk and do not consider the impact of disk migra¬ 
tion. 

The emission and absorption spectra of giant exoplan¬ 
ets and substellar companions can inform us about the 
chemical species and the C/0 ratio in their atmospheres 
and help resolve the debate about their formation mech¬ 
anism and evolution. Two complementary approaches 
have been used in constraining the atmospheric struc¬ 
ture and composition of these objects. The first one re¬ 
lies on sophisticated forward modeling of the physical 
and chemical processes in the atmospheres of giant plan¬ 
ets and calculating molecular abundances and pressure- 
temperature (P-T) profiles, and subsequently computing 
the resulting spectra (e.g., Fortney et al. 2006, 2008; Bur¬ 
rows et al. 2007, 2008). This approach, however, is com¬ 
putationally expensive, which combined with the large 
number of free parameters makes it extremely difficult 
to compare the model spectra to observations in a statis¬ 
tically robust manner. While model spectra that match 
the data can be found, quantifying the uncertainties in 
the underlying individual parameters is difficult. 

The other approach, “spectral retrieval”, on which we 
focus here, relies on a simple radiative transfer model 
that assumes that the P-T profile, the surface gravity 
and the atmospheric composition (with only few main 
molecular species considered) are free parameters (e.g., 
Benneke & Seager 2012, 2013; Lee et al. 2012; Line et al. 

2012, 2013; Heng & Lyons 2016). The effects of clouds 
on substellar emission spectra can also be parametrized 
and included (e.g., Madhusudhan et al. 2011; Lee et al. 

2013, 2014a, b). Model spectra produced in this way are 
computationally cheap and the limited number of pa¬ 
rameters makes fitting to observations feasible. Since the 
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physical phenomena included in these models are limited, 
they rely entirely on observational data to inform the fit, 
although priors based on theoretical considerations can 
improve the retrieval. A caveat is that the models used 
for retrieval, while still accounting for the most impor¬ 
tant physical processes and chemical species, could still 
be over-simplistic. Thus, not all features of a given at¬ 
mosphere may be accounted for correctly. 

In this study, we examine an integral field near-infrared 
spectrum observed by Hinkley et al. (2013) of the directly 
imaged substellar companion orbiting the young star k 
Andromedae (hereafter k And b following Carson et al. 
2013). We summarize the known properties of the k And 
system in Section 2. In Section 3, we discuss the available 
observed spectrum of k And b. We present our simple 
atmospheric model in Section 4. Section 5 focuses on our 
spectral retrieval approach, and in Section 6 we discuss 
our results and compare them to similar studies of giant 
hot exoplanets. 

2. KNOWN PROPERTIES OF THE k ANDROMEDAE 
SYSTEM 

K Andromedae is a young but already post-main se¬ 
quence B9IVn star (Wu et al. 2011) whose accurate 
age is debated. Carson et al. (2013), who discovered 
the companion, assume that k And is a member of the 
Columba association (Zuckerman et al. 2011) and, fol¬ 
lowing Marois et al. (2010), adopt an age of 30tioMyr 
for their analysis. On the other hand, Hinkley et al. 
(2013) use theoretical stellar isochrone tracks (Bertelli et 
al. 2009) combined with the previously measured stel¬ 
lar properties to derive an age of 220 ± 100 Myr. An 
alternative analysis by Bonnefoy et al. (2014) adopts 
a more conservative value than Carson et al. (2013), 
but still based on the age of the Columba association: 
30^i[g° Myr. The correct age of the system is impor¬ 
tant when assessing the luminosity and hence the mass 
of K And b, whose spectral type is LI ± 1 (Hinkley et al. 
2013). There is no available dynamical mass estimate for 
this object. Carson et al. (2013) use their estimate for 
the age of the system (30tio ^^Y^) and the DUSTY evo¬ 
lutionary models (Allard et al. 2011) to estimate a com¬ 
panion mass of 12.8'!ii'QMj. On the other hand, Hinkley 
et al. (2013) adopt Mcomp = 50 )^}|Mj, while Bonnefoy 
et al. (2014) place a lower limit on the mass of the object 
of Mcomp > lOMj based on “warm-start” evolutionary 
models. While in principle it is a basic parameter, the 
mass of the companion is of secondary importance for 
our purposes, because it enters the spectrum formation 
only via the surface gravity term, which also depends 
on the radius of the object. The effective temperature 
of K And b is 2040 ± 60K (Hinkley et al. 2013). We 
summarize the known properties of the host star and its 
companion in Table 1. 

3. OBSERVATIONS 

In our analysis, we consider the spectrum of n And b as 
published by Hinkley et al. (2013). We offer a brief sum¬ 
mary of how the data were obtained. These authors used 
the “Project 1640” instrument (Hinkley et al. 2011; Op- 
penheimer et al. 2012) on the 200 inch Hale Telescope 
at Palomar Observatory. Project 1640 is an integrated 
combination of a coronagraph and an integral field spec¬ 


trograph (IPS) and covers the Y, J and H bands. The 
data were obtained on UT 2012 December 23, starting 
at an airmass of 1.02. A total of 16 images were taken 
with exposures of 183 s each, with the Hale Telescope 
adaptive optics on, and the star hidden behind the coro- 
nagraphic mask. The spectrum extraction is described 
in detail in Hinkley et al. (2013). The observations cover 
the range between 0.9 and 1.8 ^m, but for this study, 
we focus on wavelengths longer than 1 /rm. At shorter 
wavelengths, the molecules that we consider (CH 4 , CO, 
CO 2 , H 2 O, C 2 H 2 , HCN, NH 3 , see Section 4) and the P-T 
profile have a limited impact on the emission spectrum 
(e.g., Lee et al. 2014a). The spectrum is low-resolution 
(R « 45) and contains 28 flux points in the range we 
examine. The uncertainties derived by the observers for 
every spectroscopic channel include the uncertainties due 
to photon noise, systematic errors and uncertainties re¬ 
sulting from the spectral calibration. Unfortunately, the 
absolute calibration of the resulting flux values is not pre¬ 
cise, and therefore it is difficult to derive the radius of 
the companion from its luminosity using the well known 
distance to k And (52 ± 2 pc from Hipparcos parallaxes; 
Perryman et al. 1997; Carson et al. 2013). 

4. ATMOSPHERIC EMISSION MODEL 
4.1. General Description 

For this study we require a relatively fast and there¬ 
fore relatively simple radiative transfer code that takes 
a limited number of parameters - chemical species vol¬ 
ume mixing ratios, P-T profile and a scaling factor that 
incorporates the companion’s distance and luminosity - 
and computes an emission spectrum. To this end, we im¬ 
plement a custom one-dimensional plane-parallel atmo¬ 
spheric model that we validate against the results from 
the NEMESIS code (Irwin et al. 2008) and the spectral 
synthesis software used by Line et al. (2012) and Line 
et al. (2013). These codes were designed to study the 
spectra of planets, including hot Jupiters and young di¬ 
rectly imaged gas giants. While we occasionally refer 
to them as “radiative transfer” models for simplicity, we 
do not enforce radiative equilibrium in our code, since 
K And b is a self-luminous object releasing its primor¬ 
dial heat of accretion. Due to its separation from the 
primary, the companion is negligibly externally heated. 

Following Line et al. (2012), we consider the effects of 
nine chemical species on the emission spectrum in our 
radiative transfer model: CH 4 , CO, CO 2 , H 2 O, H 2 and 
He, as well as C 2 H 2 , HCN and NH 3 . The first four 
are expected to be the major sources of molecular line 
opacity in a hot substellar atmosphere (e.g., Tinetti et 
al. 2007; Swain et al. 2009). Other relatively common 
molecules also have absorption bands in the wavelength 
range we consider, e.g. hydrogen cyanide (HCN), acety¬ 
lene (C 2 H 2 ) and ammonia (NH 3 ; e.g., Rothman et al. 
2009). The abundances of these species in self-luminous 
gas giant planets and brown dwarfs are expected to be 
much lower compared to the abundances of water, carbon 
monoxide, carbon dioxide and methane (e.g., Lodders & 
Fegley 2002; Zahnle & Marley 2014), but nevertheless 
could have an impact on the emitted spectrum and there¬ 
fore we explore models that include these three species. 
Alkali metal atomic species, e.g., K I and Na I, have some 
lines that fall within the wavelength range we consider. 
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Table 1 

Stellar and Companion Parameters for the And System 


M* (Mq) 

J* (mag)t 

H* (mag)t 

Ka,* (mag)t 

Teff,* (K) 

SpT* 

Distance (pc) 

2.4-2.5 

4.624 ± 0.264 
4.595 ±0.218 
4.571 ± 0.354 
10730 ± 250 
B9IVn 

52 ± 2 

Carson et al. (2013) 

2MASS, Skrutskie et al. (2006) 
2MASS, Skrutskie et al. (2006) 
2MASS, Skrutskie et al. (2006) 

Wu et al. (2011) 

Wu et al. (2011) 

Hipparcos] Perryman et al. (1997); 
Carson et al. (2013) 

Jcomp 

16.3 ±0.3 

Carson et al. (2013) 

Hcomp 

15.2 ±0.2 

Carson et al. (2013) 

Kg,comp 

14.6 ±0.4 

Carson et al. (2013) 

Lcomp (mag) 

13.12 ±0.09 

Carson et al. (2013) 

Teff.comp (^) 

2040 ± 60 

Hinkley et al. (2013) 

Projected separation (") 

1.064 ± 0.006 

Carson et al. (2013) 

Projected separation (AU) 

56 ± 2 

Carson et al. (2013) 

Position angle (°) 

55.9 ±0.4 

Carson et al. (2013) 


^ Two Micron All Sky Survey (2MASS) magnitude of the star (from the Infrared Science 
Archive: http://irsa.ipac.caltech.edu). 


To investigate this, we compare the k And bobserva- 
tion with an observation from the SpeX Prism Spectral 
Libraries® of a L2.5 brown dwarf with similar effective 
temperature (2MASS J14343616+2202463; Sheppard & 
Cushing 2009). In a typical example, it is possible that 
the K And b observation at 1.14/im is affected by a Na 
I line doublet. However, we would need signal to noise 
3-5 times higher to be able to make this claim. Thus, 
due to the limited resolution and large uncertainties of 
our data, neutral atomic alkali lines are unlikely to have 
a significant impact on our results. 

For CO, CO2, H2O we adopt opacities from the 
HITEMP database (Rothman et al. 2010). For methane, 
acetylene, hydrogen cyanide and ammonia, we adopt 
the opacities from Freedman et al. (2014). In addition, 
we include the collision induced opacity due to H2-H2 
and H2-He interactions (Borysow et al. 2001; Borysow 
2002). The adopted opacities cover the temperature 
range between 500 and 3000 K and the pressure range 
between 10“® and 10^ bar. The line-by-line cross section 
databases as a function of T and P are sampled at a res¬ 
olution of lcm“^. Throughout this study, we consider 
the abundance of each species as a fraction of the to¬ 
tal number of molecules. The CH4, CO, CO2 and H2O 
abundances are free parameters in all of our fits. We run 
fits with the abundances of C2H2, HCN, NH3 all fixed at 
zero or all set to be free. The ratio between H2 and He 
is fixed to 86:14 (Line et al. 2012). The mean molecu¬ 
lar weight of the atmosphere is calculated based on these 
main chemical species, assuming that the contribution 
from other molecules is small. 

In our model, the atmosphere is divided into a num¬ 
ber of horizontal layers as a function of pressure. In this 
study, we use 90 layers that are evenly spaced in log{P), 
varying between 2.4 x 10“® and 163.3 bar. Each layer 
has a temperature, based on the input P-T profile and a 
chemical composition associated with it. Following Line 
et al. (2012), we keep the relative abundances by number 
of all molecules we consider uniform for all layers, since 
the information content of the spectrum is insufficient 
to allow the retrieval of chemical gradients. We calcu¬ 
late the amount of flux from a given layer that leaves 

® http://pono.ucsd.edu/-adain/browndwarfs/spexprism/ 


the atmosphere unimpeded and integrate over all layers, 
for a given wavelength. In other words, like Line et al. 
(2012), we follow, e.g., Gray (2005) and compute, for 
each atmospheric layer, the change in optical depth, 

AP . 

= fx,k): (1) 


where AP^ is the change in pressure in the i*^*' layer (in 
units of barye, i.e., gcm“^s“^); /r (in gmol“^) and g 
(in cms^) are the mean molecular weight and the sur¬ 
face gravity, respectively. CTA.fc (in cm^) and f\^k are 
the molecular cross-section at wavelength A and mixing 
ratio of the k*^ chemical species, respectively. We inte¬ 
grate this value over all layers lying above the i*^ layer 
to calculate the total amount of optical depth, that 
a photon emitted from this layer would encounter. Since 
some of the light emitted by any given layer will be ab¬ 
sorbed within that layer, we account for this effect by 
dividing the i*** layer in upper and lower half in pressure 
by dividing APi by two and adding the Ar^^ a value of 
the upper half to ti^a- The emitted intensity of a ray of 
light through the atmosphere is then, 


h = 


90 

i=l 


P, A 

cos{9 )’ 


( 2 ) 


where B is the Planck black body radiation at tempera¬ 
ture T; and 9 is the angle between the direction to the 
observer and the normal to the surface. In order to get 
the observable flux, we integrate this over the visible sur¬ 
face of the planet, 

/*27r pTT 

Fx= I Ix cos{9)d9d(j), (3) 

Jo Jo 

To compute this integral efficiently, we take advantage 
of the Gaussian quadrature approximation. The same 
operation is performed for every wavelength in which we 
are interested between 1 and 1.8 /im. The layers at the 
top of the atmosphere have too low optical depths to have 
a significant impact on the spectrum in the near-infrared, 
while the layers near the bottom are too obscured by the 
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opacity of the upper layers (the top 10 and bottom 10 
layers have < 0 . 001 % contribution to the total flux). 

We sample the model spectrum at a resolution 50 
times higher than that of the observed spectrum - ev¬ 
ery ~ 3cm“^, or ~ OA in this wavelength range. The 
model emitted flux is binned in wavelength to match the 
resolution of the observed spectrum. This approach is 
similar to the methods of Line et al. (2012) and Line et 
al. (2013). Since the molecular opacity data is sampled 
at 1 cm“^, finer sampling than that will not improve the 
information content of our model. Sampling our model 
spectrum at full 1 cm“^ resolution slows down the model 
calculation by a factor of three compared to sampling 
at our choice of 3cm“^ and makes Bayesian fitting im¬ 
practical. To ensure robustness, we produce test model 
spectra sampled at lcm“^ and 3cm“^ for several arbi¬ 
trary but reasonable parameter combinations and find 
that, after binning the models to the resolution of the 
data, they are different by < 1 % per spectroscopic point, 
much less than the uncertainties of the observed spec¬ 
trum. Therefore, our choice of wavelength sampling of 
the models does not impact our results. 

4.2. Treatment of Clouds 

Previous emission spectrum retrieval studies have 
treated clouds and hazes (Lee et al. 2013, 2014a, b). 
These authors have explored cloudless and uniformly 
cloudy (the atmosphere being uniformly permeated by 
particles of a given size) models as well as an “interme¬ 
diate” model, where the cloud layer has a minimum and 
maximum pressure, similarly to Burrows et al. (2011). 
Lee et al. (2014a) have discussed aerosol treatment in 
detail and included it in their retrieval for the tran¬ 
sit spectrum of HD 189733b by varying the wavelength 
dependent number density of the aerosol particles and 
their size. We adopt a more numerically simplified ap¬ 
proach, choosing to parametrize the macroscopic quan¬ 
tity of cloud optical depth instead of the microscopic 
quantities of cloud particle opacity and number density. 
Since our radiative transfer code is one-dimensional, we 
consider a single cloud layer that covers the whole planet 
with no gaps. The clouds are represented by adding 
“grey opacity” to the molecular opacities in the atmo¬ 
spheric layers in the model where the cloud occurs. This 
is usually equivalent to assuming that the cloud particles 
are much larger than the observation wavelength, but 
since we consider a relatively narrow wavelength range 
between 1 and 1 . 8 /rm, there are cloud chemical species 
that would appear approximately grey even if their par¬ 
ticles are comparable in size to the observed wavelength 
(e.g., Morley et al. 2012; Lee et al. 2013). The amount 
of radiation that passes through the cloud is determined 
by the amount of extra optical depth contained in the 
“cloudy” atmospheric layers. The optical depth contri¬ 
bution from the cloud deck is invariant with wavelength. 
We describe the optical depth contribution of the cloud 
deck as a Gaussian as a function of pressure in log space, 
similarly to Heng et al. (2012). The altitude of the deck 
(i.e., the pressure at the location of the peak of the Gaus¬ 
sian), the peak optical depth contribution of the cloud 
and the vertical extent of the cloud are emission model 
inputs and can be explored by the fitting routine as free 
parameters. This scheme is illustrated in Figure 1. It 



Change in optical depth (At) 

Figure 1. We show a schematic representation of our cloud treat¬ 
ment approach for an arbitrary combination of molecular mixing 
ratios and at an arbitrary wavelength. The relation between pres¬ 
sure and optical depth is represented by the black line. The red 
line is a Gaussian describing the contribution to the change in op¬ 
tical depth due to the cloud layer. The pressure at the peak of the 
Gaussian, its amplitude and standard deviation are free parame¬ 
ters in the fit. The black points show the total change in optical 
depth in a given atmospheric layer. 

is possible that the cloud permeates a “wide” range of 
pressures, i.e., the standard deviation of the Gaussian 
optical depth contribution is so large that the Gaussian 
is in essence flat in our region of interest. This case would 
be equivalent to the uniformly cloudy model in Lee et al. 
(2013, 2014b). 

For simplicity, we do not include the effect of scat¬ 
tering that will effectively reduce the “apparent” opac¬ 
ity. Forward scattering particles will transport more light 
through the top of the cloud, which is the equivalent of 
reducing the cloud opacity in pure absorption. Thus, 
scattering and cloud opacity will be degenerate in our 
simple framework. For a detailed theoretical discussion 
of clouds in brown dwarfs and young planets, see Morley 
et al. ( 2012 ). 

5. FITTING PROCEDURE 
5.1. Free Parameters 

We explore four cases of the model fit to the observed 
spectrum - 1 ) without clouds and only including the four 
main molecules: GH 4 , GO, CO 2 and H 2 O; 2) cloudy and 
only including the four main molecules; 3) cloud-free, 
including the four main molecules plus the three addi¬ 
tional ones: C 2 H 2 , HGN, NH 3 ; and 4) cloudy and all 
seven considered molecules. We summarize the param¬ 
eters and their definition in Table 2. For each case, we 
indicate which parameters are left free and which are 
fixed at 0. The spectrum model requires that 90 layers 
in pressure and temperature are specified. However, the 
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temperature cannot be let to vary freely in each layer, 
because the spectrum does not contain sufficient infor¬ 
mation for this - even neglecting the signal-to-noise con¬ 
siderations, there would be more than 90 free parame¬ 
ters compared to a spectrum of 28 points. We select 
nine individual layers, uniformly interspersed in pressure, 
where the temperatures are free parameters (Ti — Tg). 
We then use quadratic interpolation to determine the 
temperature values for the 81 remaining layers between 
these points. In order to test whether a 9-point P-T pro¬ 
file is reasonable, we run two MCMC chains (one mil¬ 
lion iterations each) with clouds turned off and by fixing 
log{gsurf) = 4.5 using a 7-point and a 12-point P-T pro¬ 
file. We find that the molecular volume mixing ratio val¬ 
ues resulting from these are consistent with the MCMCs 
using a 9-point model well within Itr. 

The surface gravity of the companion, gsurf, is uncer¬ 
tain, and expected to be correlated with the molecular 
abundances. High surface gravity flattens the spectrum 
and requires a higher abundance of a molecular species to 
produce a feature of a given depth. Thus, hxing the sur¬ 
face gravity to a “reasonable” value would result in chem¬ 
ical abundance estimates with inaccurate uncertainties. 
Hence, even though we cannot constrain this property, 
we leave it as a free parameter in all four cases that we 
explore in detail, such that log(gsurf) G (3.5, 5.2), in cgs. 
We pick this interval based on estimates of the surface 
gravity of k, And b from Hinkley et al. (2013). 

While the relative flux from point to point within the 
spectrum is well calibrated, the absolute flux calibration 
is uncertain. To compensate for this and for the unknown 
radius of k And b, we introduce a multiplicative term, 

Acaiib ~ Gconv i where Gconv is the unknown con¬ 

version factor between measured and “true” flux. Rcomp 
is the radius of the companion and d is the distance to 
the K And system. 

5.2. Markov Chain Monte Carlo 

As a first step in the fitting process, we employ the 
IDL MPFIT y^-minimization library (Markwardt 2009) 
to provide us with an initial guess about the location of 
the minimum in the free parameter space. Then, we 
employ a Markov Chain Monte Carlo (MCMC) approach 
to get best fit values and uncertainties. For this analysis, 
we have utilized the MCMC algorithm implemented and 
described by Todorov et al. (2012, 2014), who follow Ford 
(2005, 2006). We use the Metropolis-Hastings algorithm 
within the Cibbs sampler in conjunction with a Caussian 
trial value probability distribution. We assume flat pri¬ 
ors for all parameters. During a given MCMC iteration, 
we perturb only a single randomly chosen free parameter. 
Each MCMC chain that we execute has 10® iterations, 
or approximately 50,000-70,000 iterations per parameter. 
The initial part of the MCMC chain (10% or 10® itera¬ 
tions) is discarded as “burn-in” time necessary for the 
algorithm to converge. Before running a long chain, we 
run shorter chains in order to choose the “step sizes”, 
i.e., the standard deviations of the Caussian probability 
distributions used to select the amount by which a given 
parameter is perturbed. In order to optimize the conver¬ 
gence speed of the chain, we select, by experimentation, 
the Gaussian widths such that the acceptance rate of the 
given parameter is between 15 and 35%. We constrain 


the molecular number fractions of the molecules to be 
less than 10 “^ and the temperatures are constrained to 
be less than 3000 K by preventing any trial states that 
exceed these values. This is done in order to prevent 
the chains from entering temperature regimes where the 
opacities are not defined and preventing the total abun¬ 
dance of the molecules included in the fit from becoming 
unrealistically large. 

6. RESULTS AND DISCUSSION 

6 . 1 . Adopted Parameter Values and Uncertainty 
Estimation 

As discussed in Section 5.1, we run an MCMC chain 
for each of the four cases we consider - cloudy vs. cloud- 
free and four main molecules vs. all seven species in the 
model atmosphere. For all chains, the best fit model is 
chosen as the model with the minimum y^ value in the 
MCMC chain, after the initial 10® iterations have been 
dropped as mentioned before. While this model best fits 
the available data, its parameters are not necessarily the 
values on which the MCMC converges Hnally. We find 
that the reduced y^values for the fits are 1.25, 2.21, 1.13 
and 2.06 for cases 1-4 (as dehned in Table 2), respec¬ 
tively. The corresponding Bayesian information criterion 
(BIC) values are 66.1, 82.1, 71.3 and 84.4. Even though 
the reduced y^ values of the best fit models that the 
four-species chains encountered are larger than the val¬ 
ues produced by the seven-species fits, the BIC values 
for the latter models are significantly larger. Similarly, 
cloud-free BIC values are much smaller than cloudy ones. 
This suggests that the improvements in the fit achieved 
by adding molecules and clouds are not significant given 
the extra free parameters in the model, and for our re¬ 
sults discussion, we adopt the outcome of the chain with 
the minimum BIC, case 1: cloud-free, including only 
CH 4 , CO, CO 2 and H 2 O. We compare the best-htting 
emission models from this MCMC run to the data in 
Figure 2. 

The surface gravity does not converge in any of the 
tested chains and the MCMC histograms for this pa¬ 
rameter are flat. Therefore, we are unable to place any 
constraints on it. 

We present the P-T prohles for case 1 in Figure 3, 
and for all cases in the Appendix, Figure 8 . Several 
percent of the P-T profiles tested, especially for case 3, 
cloud-free with seven molecules, exhibit sharp temper¬ 
ature inversions. We find that the affected model at¬ 
mospheres contain very little water (log(nH 2 o) < - 6 ) 
but the abundances of other molecules tend to be large 
(methane, CO and CO 2 have log(n) > —2.5)). The tem¬ 
perature inversions cause these molecules to appear in 
emission instead of absorption. Since there is no wa¬ 
ter in these models, there is a trough (in emission) near 
1.4um, while on either side the CO, CO 2 and CH 4 can 
be seen in emission. The overall flux level of the spec¬ 
trum is then compensated by the Acaiib free parameter 
(Table 2). At the resolution of the observations, this con¬ 
figuration closely resembles the observed spectrum. The 
inverted P-T prohles represent likely unphysical states of 
our vector of free parameters, considering the strength of 
the required temperature inversion, that nevertheless re¬ 
sult in spectra similar to the observed. This underscores 
the importance of high-resolution spectra of substellar 
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Table 2 

Parameters of the Fit 


Parameter 

(units) 

Case Case 2^ Case 3^ Case 4*^ 

Description 

Main molecules 


iiHaO 

nCH4 

nco 

ncOa 

Other molecules 

free 

free 

free 

free 

free 

free 

free 

free 

free 

free 

free 

free 

free 

free 

free 

free 

Water molecules number fraction. 

Methane molecules number fraction. 

CO molecules number fraction. 

CO 2 molecules number fraction. 


0 

0 

free 

free 

C 2 H 2 molecules number fraction. 

uhcn 

0 

0 

free 

free 

HCN molecules number fraction. 

iiNHa 

0 

0 

free 

free 

NH 3 molecules number fraction. 

Cloud parameters 

Camp 

0 

free 

0 

free 

Cloud peak optical depth; amplitude in Fig. 1. 

Cpress 


free 


free 

Pressure where the cloud layer is centered. 

Cvert 


free 


free 

Vertical extent of the clouds (bar). 

Physical properties 

gsurf (cms 2 ) 

free 

free 

free 

free 

Surface gravity, log(gsurf) S (3.5, 5.2). 

Ti - Tg (K) 

free 

free 

free 

free 

Temperatures throughout the atmosphere. 

Calibration 

-^calib 

free 

free 

free 

free 

Multiplicative offset (poor absolute calibration). 

^ Cloud-free, including 

only the four 

main 

molecules. 




^ Cloudy, including only the four main molecules. 
^ Cloud-free, including all seven molecules. 
Cloudy, including all seven molecules. 


companions for accurate atmospheric composition deter¬ 
mination. These states are permitted in our simplified 
model fitting framework that does not impose any re¬ 
strictions from atmospheric heat transport or chemical 
equilibrium. Were they the majority or a large fraction of 
all MCMC states, especially in case 1, we would conclude 
that the observed spectrum does not contain sufficient in¬ 
formation to constrain the atmospheric properties using 
our approach. However, they constitute a minority of the 
MCMC states. Thus, we conclude that the inverted P-T 
profile states have little influence on our final results. 

We estimate uncertainties and adopt best fit values 
based on the histograms of the MCMC parameter states. 
Even though we adopt the case 1 converged parameter 
values for our discussion based on the BIC values for 
the best ht models from the four cases, we present the 
results for both the cloudy and cloud-free models in the 
Appendix (Table 4 and Figures 8 , 9 and 10). 

We are able to estimate the volume mixing ratio of 
water molecules and we show the cloud-free water abun¬ 
dance histograms in Figure 4 for case 1, although the 
results are consistent between cases (Table 4). The 
case 1 histogram for CH 4 is flat, indicating poor con¬ 
vergence, but has a sharp peak and then cut-off to¬ 
wards higher number fractions. The cut-offs is well be¬ 
low the maximum value of 0.1 permitted in the MCMC. 
Thus, we place upper limits on the abundance of CH 4 
by reporting the number density value below which 95% 
of MCMC iterations have occurred for that parameter. 
The abundances of CO and CO 2 remain unconstrained. 
The water volume mixing ratio and the upper limit for 
methane’s volume mixing ratio are consistent with pre¬ 
dictions by Hubeny & Burrows (2007) for brown dwarfs, 
and with estimates by Heng & Lyons (2016) for C/O=0.5 


at P = 1 bar, assuming solar metallicity. 

Both CO and CO 2 have absorption bands in the 1 
to 1 . 8 /rm range, but these bands are weak, compared 
to water at these temperatures and pressures (HITEMP 
database Rothman et al. 2010), therefore, it is not sur¬ 
prising that these species, along with C 2 H 2 , HCN, NH 3 
have little impact on our results, and especially the mea¬ 
sured water volume mixing ratio. For all cloudy and 
cloud-free models, the water fractions are consistent with 
each other suggesting that clouds do not dominate the 
emission spectrum of the companion at this wavelength 
range, in the context of our simple cloud model. 

The P-T profiles we present in Figure 3 do not encom¬ 
pass all of the nine levels in pressure where the temper¬ 
ature is a free parameter. The temperatures at altitudes 
above the ~ 0.1 bar level do not converge, since these 
layers are physically located well above the part of the 
atmosphere where the optical depth is close to 1 , for all 
wavelengths between 1 and 1 . 8 /rm. 

6 . 2 . Plausibility of Grey Clouds 

Our main result, the volume mixing ratio of water is 
consistent between clear and grey cloud models, with and 
without the inclusion of additional molecules. Within 
our framework, grey clouds are not required to explain 
the observed spectrum, but cannot be excluded in prin¬ 
ciple. In order to check the physical plausibility of grey 
clouds consisting of large particles, we use typical values 
to estimate the thermal velocity of the gas in the at¬ 
mosphere and the “eddy mixing strength” (or e.g., 
Heng & Demory 2013) to be K^z > 10® — 10^®. This 
value, compared to Kz^ ^ 10"^ that Marley et al. (2012) 
find for HR 8799b, is likely unphysical, suggesting that 











The Water Abundance of k And b 


7 



Figure 2. We compare about 300 instances of the model cloud-free four-core-molecule MCMC run (orange lines) to the observed spectrum 
of K And b (connected black points; Hinkley et al. 2013). Each shown model is randomly chosen from a set with Ax^ < 1 compared to 
the model with minimum x^ = 16.1 (x^educed “ 1-25). The water absorption feature at 1.4^m is clearly seen in the spectrum. While we 
plot absolute units on the vertical axis, the absolute calibration of the spectrum was poor in the original observations (unlike the relative 
calibration at different wavelengths). Thus, the ordinate axis here should be used to compare the difference of the fluxes at different bands 
and not the absolute flux of k And b. The details of the observed spectrum and its uncertainties are discussed in Section 3. 


the large-particle, grey-cloud assumption is too simplis¬ 
tic. However, there are proposed cloud compositions for 
which the grey assumption holds approximately over the 
relatively narrow wavelength range covered by our data, 
e.g., Na 2 S (Morley et al. 2012) and MgSiOa, (Lee et al. 
2013) even for smaller particle sizes, comparable to the 
radiation wavelength. 

6.3. Comparison with Other Planets 

We compare the water abundance we retrieve for 
K And b based on case 1 with retrieved water abundances 
from other studies in Figures 6 and 7. Most of the points 
presented in these plots assume clear cloud-less models. 
Despite that, there are possible heterogeneities in the 
water volume mixing ratio determination due to differ¬ 
ent assumptions, as well as due to the different methods 
used to obtain spectral information. There are also pos¬ 
sible heterogeneities related to the temperature determi¬ 
nation. All hot Jupiters’ temperatures given in Figure 6 
are day side equilibrium temperatures based on distance 
to the star and stellar luminosity and assuming no re¬ 
distribution to the night side of the planet. For the two 
directly imaged objects, HR 8799b and k And b we adopt 
the effective temperatures found by Hinkley et al. (2013) 
and Lee et al. (2014b). 

The water volume mixing ratio in the atmosphere of 


K And b is similar to the retrieved values from the 
day side emission of the transiting hot Jupiters TrES-2b 
(1.3Mj), TrES-3b (1.9Mj) and HD 149026b (0.36Mj, 
Southworth 2010), as Eigures 6 and 7 suggest (Line et 
al. 2014). These three hot Jupiters have day side equilib¬ 
rium temperatures of ^ 1930,1920 and 2110 K, respec¬ 
tively, which is similar to the effective temperature of 
2040 ± 60K for k And b (Hinkley et al. 2013). This is 
despite the fact that k And b is likely more massive than 
all three of these, and in fact than any of the other ob¬ 
jects presented in Figures 6 and 7. The water volume 
mixing ratio of k And b falls within the range of hot 
Jupiter abundances, although the uncertainties in the 
water abundances of the transiting planets are far larger 
than those for the directly imaged planets. We attribute 
this to the fact that in most cases the model retrieval on 
secondary eclipse data relied on multiband photometry of 
the planetary day side rather than spectroscopy. In addi¬ 
tion, directly imaged spectra do not suffer from system¬ 
atic effects resulting from eclipse and transit depth mea¬ 
surements. The wide range of water abundances covered 
by the companions in Figures 6 and 7 provides no proof 
that the amount of water in a companion’s atmosphere 
is correlated with its temperature, irradiation level or 
mass. However, this may change as data with higher 
signal-to-noise ratio and with better spectral resolution 
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1000 2000 


Temperature (K) 

Figure 3. Histogram of the temperature values at each pressure 
level from the MCMC run for case 1: without clouds and only 
including CH 4 , CO, CO 2 and H 2 O. The best-fit temperatures at 
pressures corresponding to the Ti — T 9 free parameters are shown 
as red points with error bars. The uncertainties are based on the 
corresponding MCMC histograms. We define temperature param¬ 
eters with uncertainties larger than 500 K as unconstrained (not 
shown). 


Table 3 

Results from the Case 1 
MCMC Run 


Parameter^ 

Best Value'^ 

log{nH2o) 

q c-l-O.S 

log(ncH4)‘' 

log(nco) 

< - 2.6 
no constraint 

loglncos) 

no constraint 

Tg (K)d 

1430 ± 290 

Tr (K) 

1850 ± 300 

Ts (K) 

1920± 210 

e 

1.25 

BIC^ 

66.1 


^ The parameter symbols are de¬ 
fined in Table 2. 

^ Parameter values based on the 
MCMC run with no clouds and 
only the four main molecules: 
CH 4 , CO, CO 2 and H 2 O. 

^ 95% upper limits based on the 
MCMC histograms for these pa¬ 
rameters. 

^ Temperature parameters Ti 
through T 5 and Tg are uncon¬ 
strained. We define the tem¬ 
perature parameter to be “con¬ 
strained” if the 1 cr uncertain¬ 
ties from the MCMC run are less 
than 500 K. 

® The reduced value of the 
best fit model from the MCMC 
run. 

^ The Bayesian Information Cri¬ 
terion value of the best fit model 
from this MCMC run. 



Figure 4. The histogram of cloud-free MCMC parameter param¬ 
eter states for the molecular number fraction of water in log space 
for the cloud-free case, including only the four main molecules (case 
1). The MCMC histogram peaks strongly near log(nH 2 o) = —3.5. 
We denote the peak of the distribution for case 1 (solid line) and the 
region that covers 34% (“Icr”) of the values on either side (dashed 
lines). 



Figure 5. Similar to Figure 4 but for the CH 4 volume mixing ratio 
for the cloud-free, four molecules case (case 1). Despite the peak 
near log(ncH 4 ) ~ —2.5, the long tail towards very low abundances 
is pronounced, making it difficult to constrain the abundance of 
methane precisely. The red dashed line corresponds to the 95% 
upper limit on the number fraction. 

and coverage become available for the objects, allowing 
for more precise abundance measurements. 

7. SUMMARY AND FUTURE WORK 

We find that the water abundance we derive for 
K And b is not qualitatively different from these prop¬ 
erties derived from the day side spectra and photometry 
of transiting hot Jupiters, although there is a wide range 
of water abundances that are currently consistent with 
the hot Jupiter measurements. Since H 2 O is essential 
in the formation of the C/0 ratio of an atmosphere, the 
similar water abundances may point to a common forma¬ 
tion mechanism and location within the protoplanetary 
disk, followed by migration (e.g. Oberg et al. 2011; Mad- 
husudhan et al. 2014a). This suggestion can be tested as 
the quality of secondary eclipse and direct imaging ob¬ 
servations improves. Therefore, studying the molecular 
abundances of gas giant planets and companion brown 
dwarfs with very different effective temperatures, irradi¬ 
ation levels and masses that currently appear to inhabit 
a similar region of parameter space, is a logical next step. 

Future development of our computational approach 
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Planet Temperature (K) 


Figure 6. We compare the retrieved water abundances in the 
atmospheres of several hot giant planets with their temperatures. 
The retrieved water abundances of hot Jupiters based on secondary 
eclipse photometry and emission spectroscopy by Line et ah (2014) 
are shown as black circles (for WASP-43b we use the updated value 
from Kreidberg et al. 2014). The blue triangles without a label 
represent the water abundances of three planets from the Line et 
al. (2014) sample based on transit spectroscopy retrieval analysis 
by Madhusudhan et al. (2014b) and Kreidberg et al. (2015). The 
two points representing directly imaged companions (this work and 
Lee et al. 2014b) are shown as red squares. The value we show for 
K, And b is based on the cloud-free models that only include the 
abundances of CH 4 , CO, CO 2 and H 2 O as free parameters, other 
than the P-T profile and the surface gravity (case 1). Line et al. 
(2014) and Madhusudhan et al. (2014b) also make the no-clouds 
assumption, unlike Kreidberg et al. (2014) and Lee et al. (2014b). 
For this plot, we use the equilibrium dayside temperatures for the 
hot Jupiters, assuming no redistribution to their night sides, and 
the effective temperatures for the directly imaged objects (Hinkley 
et al. 2013; Lee et al. 2014b). The uncertainties of the abundance 
of water in the hot Jupiters based on secondary eclipse measure¬ 
ments are typically large since most of them have only dayside 
photometry in the mid-infrared and no spectroscopy, while the re¬ 
trievals by Madhusudhan et al. (2014b) and those for HR 8799b 
and K, And b relied on near-infrared spectra. This suggests that the 
way to maximize precision in water abundance measurements is to 
rely on spectra of directly imaged companions. The coolest “hot 
Jupiter” on the plot, GJ 436b, is in fact a Neptune-sized planet 
orbiting a red dwarf on a 2.6 day orbit. The uncertainties on the 
water abundance of HR 8799b are too small to be clearly displayed 
in this logarithmic scale. 



0.01 0.10 1.00 10.00 
Planet Mass (Mj) 


Figure 7. Similar to Figure 6 , but here we compare the water 
volume mixing ratio to the mass of the planets and substellar com¬ 
panions. We adopt the masses listed in the The Extrasolar Planets 
Encyclopaedia^, except for k And b, for which we use the mass of 
Mj found by Hinkley et al. (2013). Again, we see no obvious 
dependence of water abundance on the mass of the planet. 

^http://exoplanet.eu/ 


will include proper treatment and weighting of photo¬ 
metric data points in addition to spectroscopy, which 
will allow full usage of all available data for a given 
substellar companion. While our spectrum contains 28 
data points, which is a rich data set compared to many 
other exoplanet studies, it is still not a high-resolution, 
high-signal-to-noise spectrum. Obtaining higher qual¬ 
ity spectra of directly imaged companions using the new 
generation of high-contrast instruments like GPI (Gem¬ 
ini South, Macintosh et al. 2014) and SPHERE (VLT, 
Beuzit et al. 2008) would allow for better constraints of 
their P-T profiles in a wider range of pressures and allow 
for better quantitative comparisons to the P-T profiles 
of hot Jupiters. 

While the atmospheric emission model presented here 
is simplified and could yield unphysical results by not 
taking into account potentially important chemical 
species or physical processes, it is important to consider 
along with sophisticated, but computationally expensive, 
forward models. Our relatively fast code can account 
for and place uncertainties on the most important 
constituents of an atmosphere that shape the observed 
emission spectrum. The results from this and similar 
studies could be used as input for forward models that 
include additional physical and chemical considerations. 
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Table 4 

Results from all MCMC Runs 


Parameter^ 

Case 1*^ 

Case 2 

Case 3 


Case 4 

loglnuso) 

q C--I-0.3 

q q-|-0.5 
'^•'^-6.8 

q c:+0.4 
0.0_y_8 


q cr+O.e 

0.0_7 9 

loglncH,)'" 

log(nco) 

< - 2.6 

no constraint 

< -2.4 

no 

constraint 

no constraint 

no constraint 

no constraint 

no 

constraint 

logfncoa)'' 

no constraint 

no constraint 

< - 2.2 

no 

constraint 

loglncaH,)'" 

logfnHCNr 



< - 2.4 

< - 3.0 


< - 2.6 

< - 2.2 

logfuNHa)'" 



< - 3.0 


< - 2.6 

Ts (K)d 

1430 ± 290 

no constraint 

no constraint 

no 

constraint 

T 7 (K) 

1850 ± 300 

no constraint 

no constraint 

no 

constraint 

Tg (K) 

1920 ± 210 

no constraint 

2050 ± 230 

no 

constraint 

2 e 
''-red. 

1.25 

2.21 

1.13 


2.06 

BIC^ 

66.1 

82.1 

71.3 


84.4 


^ The parameter symbols are defined in Table 2. 

^ Cases: 1) cloud-free and only the four main molecules; 2) cloudy and only with 
the four main molecules; 3) cloud-free, including the four main molecules plus the 
three additional ones; and 4) cloudy, with all seven considered molecules. 

^ 95% upper limits based on the MCMC histograms for these parameters. 

^ Temperature parameters Ti — T 5 and T 9 are never constrained. We assume 
that a temperature parameter is “constrained” if the 1 cr uncertainties from a given 
MCMC run are less than 500 K. 

® The reduced value of the best fit model from this MCMC run. 

^ The Bayesian Information Criterion value of the best fit model from this MCMC 


1000 2000 1000 2000 1000 2000 1000 2000 
Temperature (K) Temperature (K) Temperature (K) Temperature (K) 



Figure 8 . Histograms of the temperature values at each pressure level from the MCMCs for all four considered cases: 1) without clouds 
and only including the four main molecules: CH 4 , CO, CO 2 and H 2 O; 2 ) cloudy and only including the four main molecules; 3) cloud-free, 
including the four main molecules plus the three additional ones: C 2 H 2 , HCN, NH 3 ; and 4) cloudy and all seven considered molecules. 
The red points with error bars represent the best-fit temperatures at pressures where corresponding to the Ti — Tg free parameters. The 
uncertainties are based on the corresponding MCMC histograms. We consider temperature parameters with uncertainties larger than 500 K 
be unconstrained and we do not plot these. 


APPENDIX 

MCMC RESULTS FOR ALL CASES 

Here we present the results from the MCMC runs for all four cases that we explore. Table 4 lists the best fit values 
for all converged physical parameters. Figure 8 is similar to Figure 3 but shows the temperatures for all cases, not 
just for case 1. Figures 9 and 10 show the volume mixing ratio histograms for water and methane for each of the four 
cases. Even if we had adopted case 3 or case 4 that find upper limits on the number densities of C 2 H 2 , HCN and NH 3 
these would not have been meaningful, since these molecules are not expected to be present in such large quantities in 
the atmosphere of k And b. The fact that the MCMC fits place an upper limit on these values simply confirms that 
these molecules play little role in the formation of the emission spectrum at this resolution. 
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Figure 9. The histograms of cloud-free MCMC parameter parameter states for the molecular number fraction of water in log space for 
each of the four cases we consider: cloud-free and only including the four main molecules(case 1 ; black); cloudy and only including the 
four main molecules (case 2, green); cloud-free, including all seven molecules (case 3, blue); and cloudy and including all seven molecules 
(case 4, red). All MCMC histograms peak near log(nH 2 o) ~ —3.5, but the cloudy or seven-molecule models have long tails towards very 
low abundances. We denote the peak of the distribution for case 1 (solid black line) and the region that covers 34% (“Icr”) of the values 
on either side (dashed black lines). 



log(ncHj 

Figure 10. Similar to Figure 9 but for the CH 4 volume mixing ratios. The four panels show the MCMC histograms for the four model 
cases labeled. In the cloudy cases the methane abundance is unconstrained. In the cloud-free cases, have peaks near log(nQH4) ~ —2.5, 
but their long tails towards very low abundances are pronounced. The red dashed lines correspond to the 95% upper limits on the number 
fraction of CH 4 in the cloud-free models. 
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